##### read in the sensitivity results

# set the folder names
fold_name = "D:/Research/global_refor_mac/sensitivity_runs/"
stand_name = "D:/Research/global_refor_mac/model_runs/"

# Columns (probably not all 8, but maybe): deforestation 2020-2050, reduced deforestation 2020-2050 at $50/tCO2; 
# deforestation emissions 2020-2050, reduced deforestation emissions 2020-2050 at $50/tCO2; 
# reforestation 2020-2050, enhanced reforestation 2020-2050 at $50/tCO2; 
# reforestation removals 2020-2050, enhanced reforestation removals 2020-2050 at $50/tCO2.


model_calc <- function(model_list, name_vec){
 

  # deforestation 
  def_skm_2020_2050_bau = (sum(apply(model_list[[1]][[2]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000
  def_emiss_2020_2050_bau = sum(apply(model_list[[1]][[3]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000
  red_def_skm_2020_2050_10dol = def_skm_2020_2050_bau- ((sum(apply(model_list[[3]][[2]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000)
  red_emiss_2020_2050_10dol = def_emiss_2020_2050_bau -(sum(apply(model_list[[3]][[3]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000)
  red_def_skm_2020_2050_20dol = def_skm_2020_2050_bau- ((sum(apply(model_list[[4]][[2]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000)
  red_emiss_2020_2050_20dol = def_emiss_2020_2050_bau -(sum(apply(model_list[[4]][[3]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000)
  red_def_skm_2020_2050_50dol = def_skm_2020_2050_bau- ((sum(apply(model_list[[7]][[2]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000)
  red_emiss_2020_2050_50dol = def_emiss_2020_2050_bau -(sum(apply(model_list[[7]][[3]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000)
  red_def_skm_2020_2050_100dol = def_skm_2020_2050_bau- ((sum(apply(model_list[[12]][[2]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000)
  red_emiss_2020_2050_100dol = def_emiss_2020_2050_bau -(sum(apply(model_list[[12]][[3]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000)
  
  # reforestation 
  refor_skm_2020_2050_bau = (sum(apply(model_list[[1]][[5]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000
  refor_emiss_2020_2050_bau = sum(apply(model_list[[1]][[6]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000
  red_refor_skm_2020_2050_10dol = ((sum(apply(model_list[[3]][[5]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000) - refor_skm_2020_2050_bau
  red_refor_emiss_2020_2050_10dol = (sum(apply(model_list[[3]][[6]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000) - refor_emiss_2020_2050_bau
  red_refor_skm_2020_2050_20dol = ((sum(apply(model_list[[4]][[5]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000) - refor_skm_2020_2050_bau
  red_refor_emiss_2020_2050_20dol = (sum(apply(model_list[[4]][[6]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000) - refor_emiss_2020_2050_bau
  red_refor_skm_2020_2050_50dol = ((sum(apply(model_list[[7]][[5]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000) - refor_skm_2020_2050_bau
  red_refor_emiss_2020_2050_50dol = (sum(apply(model_list[[7]][[6]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000) - refor_emiss_2020_2050_bau
  red_refor_skm_2020_2050_100dol = ((sum(apply(model_list[[12]][[5]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000) - refor_skm_2020_2050_bau
  red_refor_emiss_2020_2050_100dol = (sum(apply(model_list[[12]][[6]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000) - refor_emiss_2020_2050_bau
  
  # create data frame for 10 dollars
  dollar_10_df = data.frame("name" = name_vec,
                            "Deforestation 2020-2050 (million ha)" = def_skm_2020_2050_bau,
                            "Emissions from deforestation 2020-2050 (tCO2)" = def_emiss_2020_2050_bau,
                            "Redudced deforestation at $10 (million ha)" = red_def_skm_2020_2050_10dol,
                            "Reduced emissions from deforestation at $10 (tCO2)" = red_emiss_2020_2050_10dol,
                            "Reforestation 2020-2050 (million ha)" = refor_skm_2020_2050_bau,
                            "Removal from reforestation 2020-2050 (tCO2)" = refor_emiss_2020_2050_bau,
                            "Enhanced reforestation at $10 (million ha)" = red_refor_skm_2020_2050_10dol,
                            "Enhanced removals from reforestation at $10 (tCO2)" = red_refor_emiss_2020_2050_10dol,
                            stringsAsFactors = FALSE)
  colnames(dollar_10_df) = c("name", "Deforestation 2020-2050 (million ha)", 
                            "Emissions from deforestation 2020-2050 (tCO2)", 
                             "Redudced deforestation at $10 (million ha)", 
                              "Reduced emissions from deforestation at $10 (tCO2)", 
                              "Reforestation 2020-2050 (million ha)", 
                              "Removal from reforestation 2020-2050 (tCO2)", 
                              "Enhanced reforestation at $10 (million ha)", 
                              "Enhanced removals from reforestation at $10 (tCO2)")
  
  
  # create data frame for 20 dollars
  dollar_20_df = data.frame("name" = name_vec,
                            "Deforestation 2020-2050 (million ha)" = def_skm_2020_2050_bau,
                            "Emissions from deforestation 2020-2050 (tCO2)" = def_emiss_2020_2050_bau,
                            "Redudced deforestation at $20 (million ha)" = red_def_skm_2020_2050_20dol,
                            "Reduced emissions from deforestation at $20 (tCO2)" = red_emiss_2020_2050_20dol,
                            "Reforestation 2020-2050 (million ha)" = refor_skm_2020_2050_bau,
                            "Removal from reforestation 2020-2050 (tCO2)" = refor_emiss_2020_2050_bau,
                            "Enhanced reforestation at $20 (million ha)" = red_refor_skm_2020_2050_20dol,
                            "Enhanced removals from reforestation at $20 (tCO2)" = red_refor_emiss_2020_2050_20dol,
                            stringsAsFactors = FALSE)
  colnames(dollar_20_df) = c("name", "Deforestation 2020-2050 (million ha)", 
                             "Emissions from deforestation 2020-2050 (tCO2)", 
                             "Redudced deforestation at $20 (million ha)", 
                             "Reduced emissions from deforestation at $20 (tCO2)", 
                             "Reforestation 2020-2050 (million ha)", 
                             "Removal from reforestation 2020-2050 (tCO2)", 
                             "Enhanced reforestation at $20 (million ha)", 
                             "Enhanced removals from reforestation at $20 (tCO2)")
  
  # create data frame for 50 dollars
  dollar_50_df = data.frame("name" = name_vec,
                            "Deforestation 2020-2050 (million ha)" = def_skm_2020_2050_bau,
                            "Emissions from deforestation 2020-2050 (tCO2)" = def_emiss_2020_2050_bau,
                            "Redudced deforestation at $50 (million ha)" = red_def_skm_2020_2050_50dol,
                            "Reduced emissions from deforestation at $50 (tCO2)" = red_emiss_2020_2050_50dol,
                            "Reforestation 2020-2050 (million ha)" = refor_skm_2020_2050_bau,
                            "Removal from reforestation 2020-2050 (tCO2)" = refor_emiss_2020_2050_bau,
                            "Enhanced reforestation at $50 (million ha)" = red_refor_skm_2020_2050_50dol,
                            "Enhanced removals from reforestation at $50 (tCO2)" = red_refor_emiss_2020_2050_50dol,
                            stringsAsFactors = FALSE)
  colnames(dollar_50_df) = c("name", "Deforestation 2020-2050 (million ha)", 
                             "Emissions from deforestation 2020-2050 (tCO2)", 
                             "Redudced deforestation at $50 (million ha)", 
                             "Reduced emissions from deforestation at $50 (tCO2)", 
                             "Reforestation 2020-2050 (million ha)", 
                             "Removal from reforestation 2020-2050 (tCO2)", 
                             "Enhanced reforestation at $50 (million ha)", 
                             "Enhanced removals from reforestation at $50 (tCO2)")
  
  
  # create data frame for 100 dollars
  dollar_100_df = data.frame("name" = name_vec,
                            "Deforestation 2020-2050 (million ha)" = def_skm_2020_2050_bau,
                            "Emissions from deforestation 2020-2050 (tCO2)" = def_emiss_2020_2050_bau,
                            "Redudced deforestation at $100 (million ha)" = red_def_skm_2020_2050_100dol,
                            "Reduced emissions from deforestation at $100 (tCO2)" = red_emiss_2020_2050_100dol,
                            "Reforestation 2020-2050 (million ha)" = refor_skm_2020_2050_bau,
                            "Removal from reforestation 2020-2050 (tCO2)" = refor_emiss_2020_2050_bau,
                            "Enhanced reforestation at $100 (million ha)" = red_refor_skm_2020_2050_100dol,
                            "Enhanced removals from reforestation at $100 (tCO2)" = red_refor_emiss_2020_2050_100dol,
                            stringsAsFactors = FALSE)
  colnames(dollar_100_df) = c("name", "Deforestation 2020-2050 (million ha)", 
                             "Emissions from deforestation 2020-2050 (tCO2)", 
                             "Redudced deforestation at $100 (million ha)", 
                             "Reduced emissions from deforestation at $100 (tCO2)", 
                             "Reforestation 2020-2050 (million ha)", 
                             "Removal from reforestation 2020-2050 (tCO2)", 
                             "Enhanced reforestation at $100 (million ha)", 
                             "Enhanced removals from reforestation at $100 (tCO2)")
  
  list_out = list(dollar_10_df, dollar_20_df, dollar_50_df, dollar_100_df)
  return(list_out)
  
}


model_calc2 <- function(model_list, name_vec){

  
  # deforestation 
  def_skm_2020_2050_bau = (sum(apply(model_list[[1]][[2]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000
  def_emiss_2020_2050_bau = sum(apply(model_list[[1]][[3]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000
  red_def_skm_2020_2050_10dol = def_skm_2020_2050_bau- ((sum(apply(model_list[[3]][[2]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000)
  red_emiss_2020_2050_10dol = def_emiss_2020_2050_bau -(sum(apply(model_list[[3]][[3]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000)
  red_def_skm_2020_2050_20dol = def_skm_2020_2050_bau- ((sum(apply(model_list[[4]][[2]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000)
  red_emiss_2020_2050_20dol = def_emiss_2020_2050_bau -(sum(apply(model_list[[4]][[3]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000)
  red_def_skm_2020_2050_50dol = def_skm_2020_2050_bau- ((sum(apply(model_list[[7]][[2]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000)
  red_emiss_2020_2050_50dol = def_emiss_2020_2050_bau -(sum(apply(model_list[[7]][[3]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000)
  red_def_skm_2020_2050_100dol = def_skm_2020_2050_bau- ((sum(apply(model_list[[12]][[2]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000)
  red_emiss_2020_2050_100dol = def_emiss_2020_2050_bau -(sum(apply(model_list[[12]][[3]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000)
  
  # reforestation 
  refor_skm_2020_2050_bau = (sum(apply(model_list[[1]][[4]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000
  refor_emiss_2020_2050_bau = sum(apply(model_list[[1]][[5]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000
  red_refor_skm_2020_2050_10dol = ((sum(apply(model_list[[3]][[4]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000) - refor_skm_2020_2050_bau
  red_refor_emiss_2020_2050_10dol = (sum(apply(model_list[[3]][[5]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000) - refor_emiss_2020_2050_bau
  red_refor_skm_2020_2050_20dol = ((sum(apply(model_list[[4]][[4]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000) - refor_skm_2020_2050_bau
  red_refor_emiss_2020_2050_20dol = (sum(apply(model_list[[4]][[5]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000) - refor_emiss_2020_2050_bau
  red_refor_skm_2020_2050_50dol = ((sum(apply(model_list[[7]][[4]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000) - refor_skm_2020_2050_bau
  red_refor_emiss_2020_2050_50dol = (sum(apply(model_list[[7]][[5]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000) - refor_emiss_2020_2050_bau
  red_refor_skm_2020_2050_100dol = ((sum(apply(model_list[[12]][[4]][,2:4], 2, sum, na.rm = TRUE)) * 100) / 1000000) - refor_skm_2020_2050_bau
  red_refor_emiss_2020_2050_100dol = (sum(apply(model_list[[12]][[5]][,2:4], 2, sum, na.rm = TRUE)) / 1000000000) - refor_emiss_2020_2050_bau
  
  # create data frame for 10 dollars
  dollar_10_df = data.frame("name" = name_vec,
                            "Deforestation 2020-2050 (million ha)" = def_skm_2020_2050_bau,
                            "Emissions from deforestation 2020-2050 (tCO2)" = def_emiss_2020_2050_bau,
                            "Redudced deforestation at $10 (million ha)" = red_def_skm_2020_2050_10dol,
                            "Reduced emissions from deforestation at $10 (tCO2)" = red_emiss_2020_2050_10dol,
                            "Reforestation 2020-2050 (million ha)" = refor_skm_2020_2050_bau,
                            "Removal from reforestation 2020-2050 (tCO2)" = refor_emiss_2020_2050_bau,
                            "Enhanced reforestation at $10 (million ha)" = red_refor_skm_2020_2050_10dol,
                            "Enhanced removals from reforestation at $10 (tCO2)" = red_refor_emiss_2020_2050_10dol,
                            stringsAsFactors = FALSE)
  colnames(dollar_10_df) = c("name", "Deforestation 2020-2050 (million ha)", 
                             "Emissions from deforestation 2020-2050 (tCO2)", 
                             "Redudced deforestation at $10 (million ha)", 
                             "Reduced emissions from deforestation at $10 (tCO2)", 
                             "Reforestation 2020-2050 (million ha)", 
                             "Removal from reforestation 2020-2050 (tCO2)", 
                             "Enhanced reforestation at $10 (million ha)", 
                             "Enhanced removals from reforestation at $10 (tCO2)")
  
  
  # create data frame for 20 dollars
  dollar_20_df = data.frame("name" = name_vec,
                            "Deforestation 2020-2050 (million ha)" = def_skm_2020_2050_bau,
                            "Emissions from deforestation 2020-2050 (tCO2)" = def_emiss_2020_2050_bau,
                            "Redudced deforestation at $20 (million ha)" = red_def_skm_2020_2050_20dol,
                            "Reduced emissions from deforestation at $20 (tCO2)" = red_emiss_2020_2050_20dol,
                            "Reforestation 2020-2050 (million ha)" = refor_skm_2020_2050_bau,
                            "Removal from reforestation 2020-2050 (tCO2)" = refor_emiss_2020_2050_bau,
                            "Enhanced reforestation at $20 (million ha)" = red_refor_skm_2020_2050_20dol,
                            "Enhanced removals from reforestation at $20 (tCO2)" = red_refor_emiss_2020_2050_20dol,
                            stringsAsFactors = FALSE)
  colnames(dollar_20_df) = c("name", "Deforestation 2020-2050 (million ha)", 
                             "Emissions from deforestation 2020-2050 (tCO2)", 
                             "Redudced deforestation at $20 (million ha)", 
                             "Reduced emissions from deforestation at $20 (tCO2)", 
                             "Reforestation 2020-2050 (million ha)", 
                             "Removal from reforestation 2020-2050 (tCO2)", 
                             "Enhanced reforestation at $20 (million ha)", 
                             "Enhanced removals from reforestation at $20 (tCO2)")
  
  # create data frame for 50 dollars
  dollar_50_df = data.frame("name" = name_vec,
                            "Deforestation 2020-2050 (million ha)" = def_skm_2020_2050_bau,
                            "Emissions from deforestation 2020-2050 (tCO2)" = def_emiss_2020_2050_bau,
                            "Redudced deforestation at $50 (million ha)" = red_def_skm_2020_2050_50dol,
                            "Reduced emissions from deforestation at $50 (tCO2)" = red_emiss_2020_2050_50dol,
                            "Reforestation 2020-2050 (million ha)" = refor_skm_2020_2050_bau,
                            "Removal from reforestation 2020-2050 (tCO2)" = refor_emiss_2020_2050_bau,
                            "Enhanced reforestation at $50 (million ha)" = red_refor_skm_2020_2050_50dol,
                            "Enhanced removals from reforestation at $50 (tCO2)" = red_refor_emiss_2020_2050_50dol,
                            stringsAsFactors = FALSE)
  colnames(dollar_50_df) = c("name", "Deforestation 2020-2050 (million ha)", 
                             "Emissions from deforestation 2020-2050 (tCO2)", 
                             "Redudced deforestation at $50 (million ha)", 
                             "Reduced emissions from deforestation at $50 (tCO2)", 
                             "Reforestation 2020-2050 (million ha)", 
                             "Removal from reforestation 2020-2050 (tCO2)", 
                             "Enhanced reforestation at $50 (million ha)", 
                             "Enhanced removals from reforestation at $50 (tCO2)")
  
  
  # create data frame for 100 dollars
  dollar_100_df = data.frame("name" = name_vec,
                             "Deforestation 2020-2050 (million ha)" = def_skm_2020_2050_bau,
                             "Emissions from deforestation 2020-2050 (tCO2)" = def_emiss_2020_2050_bau,
                             "Redudced deforestation at $100 (million ha)" = red_def_skm_2020_2050_100dol,
                             "Reduced emissions from deforestation at $100 (tCO2)" = red_emiss_2020_2050_100dol,
                             "Reforestation 2020-2050 (million ha)" = refor_skm_2020_2050_bau,
                             "Removal from reforestation 2020-2050 (tCO2)" = refor_emiss_2020_2050_bau,
                             "Enhanced reforestation at $100 (million ha)" = red_refor_skm_2020_2050_100dol,
                             "Enhanced removals from reforestation at $100 (tCO2)" = red_refor_emiss_2020_2050_100dol,
                             stringsAsFactors = FALSE)
  colnames(dollar_100_df) = c("name", "Deforestation 2020-2050 (million ha)", 
                              "Emissions from deforestation 2020-2050 (tCO2)", 
                              "Redudced deforestation at $100 (million ha)", 
                              "Reduced emissions from deforestation at $100 (tCO2)", 
                              "Reforestation 2020-2050 (million ha)", 
                              "Removal from reforestation 2020-2050 (tCO2)", 
                              "Enhanced reforestation at $100 (million ha)", 
                              "Enhanced removals from reforestation at $100 (tCO2)")
  
  list_out = list(dollar_10_df, dollar_20_df, dollar_50_df, dollar_100_df)
  return(list_out)
  
}


# read in the standard 
standard_df = model_calc(model_list = readRDS(paste0(stand_name, "glob_model_run_stand.rds")),
                         name_vec = "Base Scenario")
gc()

# read in 10 percent forest cover
percent_for_10 = model_calc(model_list = readRDS(paste0(fold_name, "forest_cover_10percent.rds")),
                            name_vec = "10 percent forest cover")

# read in 50 percent forest cover
percent_for_50 = model_calc(model_list = readRDS(paste0(fold_name, "forest_cover_50percent.rds")),
                         name_vec = "50 percent forest cover")

# read in the 3 percent discount rate
discount_3per = model_calc(model_list = readRDS(paste0(fold_name, "discount_rate_3per.rds")),
                           name_vec = "Discount rate (3 percent)")

discount_7per = model_calc(model_list = readRDS(paste0(fold_name, "discount_rate_7per.rds")),
                           name_vec = "Discount rate (7 percent)")

# years of carbon repayments
carbon_payments_50year = model_calc(model_list = readRDS(paste0(fold_name, "years_50_carbon_repayments.rds")),
                                    name_vec = "50 years of carbon repayments")

carbon_payments_100year = model_calc(model_list = readRDS(paste0(fold_name, "years_100_carbon_repayments.rds")),
                                    name_vec = "100 years of carbon repayments")

# future agricultural prices 
future_ag_prices_05x = model_calc(model_list = readRDS(paste0(fold_name, "future_prices_x05.rds")),
                                  name_vec = "Future agricultural prices x 0.5")

future_ag_prices_2x = model_calc(model_list = readRDS(paste0(fold_name, "future_prices_x2.rds")),
                                  name_vec = "Future agricultural prices x 2")

# ag_rev_sensx05
ag_price_salience_05x = model_calc(model_list = readRDS(paste0(fold_name, "ag_rev_sensx05.rds")),
                               name_vec = "Agricultural prices salience x 0.5")

ag_price_salience_2x = model_calc(model_list = readRDS(paste0(fold_name, "ag_rev_sensx2.rds")),
                                   name_vec = "Agricultural prices salience x 2")

# median natural forest regrowth
median_natural_forests = model_calc(model_list = readRDS(paste0(fold_name, "median_age_of_plantation.rds")),
                                   name_vec = "Plantation vs natural (median value)")

all_plantation = model_calc(model_list = readRDS(paste0(fold_name, "all_plantation_percentage.rds")),
                                    name_vec = "All plantation regrowth")

all_natural = model_calc(model_list = readRDS(paste0(fold_name, "all_natural_percentage.rds")),
                         name_vec = "All natural regrowth")

# plantation rotation age 
plantation_rotation_age20 = model_calc(model_list = readRDS(paste0(fold_name, "plant_rotation_20year.rds")),
                         name_vec = "Plantation 20 year rotation age")

# global leakage parameter 
leakage_parameter = model_calc(model_list = readRDS(paste0(fold_name, "refor_leakage_model.rds")),
                                       name_vec = "Leakage of reforestation")

# global leakage parameter 
leakage_parameter_v2 = model_calc(model_list = readRDS(paste0(fold_name, "refor_leakage_model_optionb.rds")),
                               name_vec = "Leakage of reforestation")


# Net transaction costs, management costs, and co-benefits
negative_1000d = model_calc(model_list = readRDS(paste0(fold_name, "negative_hurdle_1000d.rds")),
                            name_vec = "Management costs and co-benefits (-$1000 / hec)")

positive_1000d = model_calc(model_list = readRDS(paste0(fold_name, "positive_hurdle_1000d.rds")),
                            name_vec = "Management costs and co-benefits (+$1000 / hec)")


maximum_overlap = model_calc2(model_list = readRDS(paste0(fold_name, "maximum_overlap_model.rds")),
                            name_vec = "Maximum overlap in grid cell")

oil_palm_removed = model_calc(model_list = readRDS(paste0(fold_name, "palm_oil_exclusion.rds")),
                             name_vec = "Exclude palm oil forest cover change")

hardwood_export_prices = model_calc(model_list = readRDS(paste0(fold_name, "hardwood_regression.rds")),
                              name_vec = "Include hardwood export prices")

cont_specific_model = model_calc(model_list = readRDS(paste0(fold_name, "glob_model_run_by_cont.rds")),
                                    name_vec = "Continent specific model")

count_biome_fixed = model_calc(model_list = readRDS(paste0(fold_name, "country_ecoregion_fixed_effect.rds")),
                                 name_vec = "Country and biome fixed effect model")


### calculate tweedie, negative binom, 15 perc, double 
tweedie_regression = model_calc(model_list = readRDS(paste0(stand_name, "glob_model_run_tweedie.rds")),
                               name_vec = "Tweedie specification model")

negbin_regression = model_calc(model_list = readRDS(paste0(stand_name, "glob_model_run_negbin.rds")),
                                name_vec = "Negative binomial specification model")


discount_15per = model_calc(model_list = readRDS(paste0(fold_name, "discount_rate_15per.rds")),
                           name_vec = "Discount rate (15 percent)")


simult_model = model_calc(model_list = readRDS(paste0(stand_name, "glob_model_run_simul_carbon.rds")),
                          name_vec = "Simultenous model")


require(dplyr)
require(tidyr)

# create 10 dollar spreadsheet
dollar_10d = bind_rows(standard_df[[1]], percent_for_10[[1]], percent_for_50[[1]],
          discount_3per[[1]], discount_7per[[1]], carbon_payments_50year[[1]],
          carbon_payments_100year[[1]], future_ag_prices_05x[[1]], future_ag_prices_2x[[1]],
          ag_price_salience_05x[[1]], ag_price_salience_2x[[1]], median_natural_forests[[1]],
          all_plantation[[1]], all_natural[[1]], plantation_rotation_age20[[1]],
          leakage_parameter[[1]], negative_1000d[[1]], positive_1000d[[1]],
          maximum_overlap[[1]], oil_palm_removed[[1]], hardwood_export_prices[[1]],
          cont_specific_model[[1]], count_biome_fixed[[1]], tweedie_regression[[1]],
          negbin_regression[[1]], discount_15per[[1]], simult_model[[1]])

dollar_20d = bind_rows(standard_df[[2]], percent_for_10[[2]], percent_for_50[[2]],
                       discount_3per[[2]], discount_7per[[2]], carbon_payments_50year[[2]],
                       carbon_payments_100year[[2]], future_ag_prices_05x[[2]], future_ag_prices_2x[[2]],
                       ag_price_salience_05x[[2]], ag_price_salience_2x[[2]], median_natural_forests[[2]],
                       all_plantation[[2]], all_natural[[2]], plantation_rotation_age20[[2]],
                       leakage_parameter[[2]], negative_1000d[[2]], positive_1000d[[2]],
                       maximum_overlap[[2]], oil_palm_removed[[2]], hardwood_export_prices[[2]],
                       cont_specific_model[[2]], count_biome_fixed[[2]], tweedie_regression[[2]],
                       negbin_regression[[2]], discount_15per[[2]], simult_model[[2]])

dollar_50d = bind_rows(standard_df[[3]], percent_for_10[[3]], percent_for_50[[3]],
                       discount_3per[[3]], discount_7per[[3]], carbon_payments_50year[[3]],
                       carbon_payments_100year[[3]], future_ag_prices_05x[[3]], future_ag_prices_2x[[3]],
                       ag_price_salience_05x[[3]], ag_price_salience_2x[[3]], median_natural_forests[[3]],
                       all_plantation[[3]], all_natural[[3]], plantation_rotation_age20[[3]],
                       leakage_parameter[[3]], negative_1000d[[3]], positive_1000d[[3]],
                       maximum_overlap[[3]], oil_palm_removed[[3]], hardwood_export_prices[[3]],
                       cont_specific_model[[3]], count_biome_fixed[[3]], tweedie_regression[[3]],
                       negbin_regression[[3]], discount_15per[[3]], simult_model[[3]])


dollar_100d = bind_rows(standard_df[[4]], percent_for_10[[4]], percent_for_50[[4]],
                       discount_3per[[4]], discount_7per[[4]], carbon_payments_50year[[4]],
                       carbon_payments_100year[[4]], future_ag_prices_05x[[4]], future_ag_prices_2x[[4]],
                       ag_price_salience_05x[[4]], ag_price_salience_2x[[4]], median_natural_forests[[4]],
                       all_plantation[[4]], all_natural[[4]], plantation_rotation_age20[[4]],
                       leakage_parameter[[4]], negative_1000d[[4]], positive_1000d[[4]],
                       maximum_overlap[[4]], oil_palm_removed[[4]], hardwood_export_prices[[4]],
                       cont_specific_model[[4]], count_biome_fixed[[4]], tweedie_regression[[4]],
                       negbin_regression[[4]], discount_15per[[4]], simult_model[[4]])

write.csv(dollar_10d, "D:/Research/global_refor_mac/graphics/sensitivity_analysis_10dollars.csv")
write.csv(dollar_20d, "D:/Research/global_refor_mac/graphics/sensitivity_analysis_20dollars.csv")
write.csv(dollar_50d, "D:/Research/global_refor_mac/graphics/sensitivity_analysis_50dollars.csv")
write.csv(dollar_100d, "D:/Research/global_refor_mac/graphics/sensitivity_analysis_100dollars.csv")

d1 = data.frame("category" = dollar_10d$name,
  "Reforestation 2020-2050 (million ha)"  = dollar_10d$`Reforestation 2020-2050 (million ha)`, 
"Removal from reforestation 2020-2050 (tCO2)"  = dollar_10d$`Removal from reforestation 2020-2050 (tCO2)`, 
"Enhanced removals from reforestation at $20 (tCO2)"  = dollar_20d$`Enhanced removals from reforestation at $20 (tCO2)`, 
"Enhanced removals from reforestation at $50 (tCO2)"  = dollar_50d$`Enhanced removals from reforestation at $50 (tCO2)`,
"Reduced emissions from deforestation at $20 (tCO2)" = dollar_20d$`Reduced emissions from deforestation at $20 (tCO2)`, 
stringsAsFactors = FALSE)

write.csv(d1, "D:/Research/global_refor_mac/graphics/sensitivity_print_v3.csv")


